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Abstract 

The radio frequency spectrum of the fermions in the unitary limit at finite temperatures is char- 
acterized by the sum rule relations. We consider a simple picture where the atoms are removed by 
radio frequency excitations from the strongly interacting states into a state of negligible interaction. 
We calculate the moments of the response function in the range of temperature 0.08ei? < T < 0.8ei? 
using auxiliary field Monte Carlo technique (AFMC) in which continuum auxiliary fields with a 
density dependent shift are used. We estimate the effects of superfluid pairing from the clock shift. 
We find a qualitative agreement with the pairing gap - pseudogap transition behavior. We also 
find within the quasiparticle picture that in order for the gap to come into quantitative agreement 
with the previously known value at T = , the effective mass has to be m* ~ 1.43m. Finally, we 
discuss implications for the adiabatic sweep of the resonant magnetic field. 

PACS numbers: 03.75. Ss, 05.30.Fk, 21.65.-f, 87.55.kh 
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Introduction: Experiments with the dilute fermionic atomic gases such as those of Li 6 and 
K 40 atoms have seen great developments l| . Dilute Fermi gases provide a clean and control- 
lable model for understanding different problems in a wide range of many body physics. Two 
species Fermi gas in the unitarity regime sits right in the middle of the BCS-BEC crossover 
regime and has been the focus of great interest. Radio frequency (RF) experiments can 
provide information on the interaction of the fermions by inducing excitations in the Fermi 



ew years, experimen- 
2J|. However, in these 



gas: in particular, the effects of fermionic superfluidity. Over the last 
tal measurements without the line broadening have become possible 
early experiments the final state interaction effects seem non-negligible and they are rather 
hard to interpret for T > T c : in the experiment by Chin et al. [2J, only the sharp free atom 
transitions are observed. Recently, experimental data with negligible final state interactions 
and small collisional effects have been made available (see Ref. for experiment and 

Ref. jlj for theory). In particular, the double peak structure of the RF response J4] of the 
polarized Fermi gas is qualitatively consistent with the phase separation. In this article, we 
study the RF response spectrum and the effective mass relevant for these experiments by 
using an ab initio method. Here, we consider a situation where the atoms in the excited state 
remain non-interacting (this is the likely picture for K 40 around the Feshbach resonance at 
202 G jsl). This simplified situation allows us to calculate the clock shift and to estimate 
the effects of the energy excitation gap across the T c . 

One of the most puzzling aspects of the strongly interacting fermions is the temperature 
dependence of the quasiparticle gap A qp . Previous mean field analysis as well as the recent 
ab initio calculations jf], B, Is, 9] suggest the existence of the so-called pseudogap above the 



T c (of which value at unitarity is most likely ~ 0.15ei? 10, lU[, £f = Fermi energy). The 



gap A qp in the quasiparticle spectrum as a function of the T was shown to have a slight dip 
around T ~ T c p, |9J and then to stay non-zero even at T > T c (at least up to T ~ 6f)- This 
behavior is in sharp contrast to the condensate fraction calculated by two-body correlation 
Aq = lim y(J dri J dr 2 ^\(ri + r)^\{r2 + ^)^2(^2)^i( r i)) that becomes non-zero only for 
T < T c . Although the interpretation of A qp as the order parameter does not apply for all 
temperatures, A qp is to show the transition from the superfluid pairing gap to the pseudogap 
or insulator gap for T > T c . We present in this article the temperature dependence of the RF 
spectrum in the unitarity regime in a simplified picture where the final state interactions can 
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be ignored and the states remain sharply defined. We test the accuracy of the quasiparticle 
gap measurement from the RF response spectrum. Finally, we relate the RF response to the 



adiabatic effect 



12J at T f» T c and discuss its experimental implications 



'or the observation 



of temperature changes during adiabatic sweep of the resonant field 
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Model: The scenario we consider is the one in which the trapped atoms can occupy the 
Zeeman states labeled by a — 1, 2, 3. Initially, only the states 1 and 2 are occupied by the 
atoms in equal number. The atoms in the states 1 and 2 are coupled with interaction given 
by the s-wave scattering length a s . Our model is implemented in a cubic lattice box of 
coordinate points spaced by 81 in each direction. Thus, a cutoff in the momentum k c = tt/SI 



is imposed and the coup 
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ing strength is regulated by the relation: l/g = m/(AnTi a s ) — 
111 ] where = h 2 k 2 /(2m) and f2 = volume of the system. In the 



i/nEfe =0 i/(2c k ) 

absence of the RF perturbation, the state 3 is not coupled to any of the other two states 
(unlike in the Ref. [3]). This model is captured by the linear response theory of the 
unperturbed thermodynamic potential 



n 



3 

^ (7=1 



h 2 v 2 

2m 



+ / drg^\(r)^i(r)^l(r)^f 2 (r) 



3 

<T=1 



(1) 



with the time dependent external perturbation TL pert (t). Here, the ^J-(r) and ^o-(r) are the 
usual Fermi operators. \i a are the chemical potentials corresponding to the Zeeman levels 
with energies e a (e 3 > e 2 > ei). When the state 2 is coupled to the state 3 by the RF 
perturbation, we assume a time dependent perturbation of the form TC per t(t) = y cos {uit) e 01 
with 7] being a small and positive number and y = % J (irl^I/^r) 1 !^ 1 ") — ^2( r )^3( r )] is 
a hermitian operator. The response function at any T with respect to the perturbation 
operator y is given by 

1 



S(u>) 



-E k f3 



J- e - EmP (m\y*\n)(n\y\m)5(u - u n 



(2) 



n,m 



where \n) and E n are the eigenstates and eigenvalues of the unperturbed Hamiltonian (Eq. 
[TJ). And, we have f3 = 1/T and uo nm = (E n — E m )/H. Although methods to calculate directly 
the response function at zero temperature exist in the context of the Hartree-Fock theory 



and its density functional extension [if]], to obtain the response function in the full energy 
spectrum by an ab initio method remains as an open problem. Instead, we calculate the 
moments of the response function characterizing the energy transfer by the RF signal. The 
moment of i-th power is defined as Af< = h i+1 J duu^u). Mi can be evaluated by direct 
application of Eq. [2J The sum rules in compact notation are given as M Q = Tr[3^ye~ W/3 ]/Z, 
M l = Tr[yi[H,y]e- n P]/Z, and M 2 = Tr[[y t ,7Y] [Tt, y]e~ ni3 ]/Z. Here, we have defined the 
partition function Z = Tr[e~ W/9 ]. These thermal averages are evaluated in the unperturbed 
basis where the system has unoccupied state 3. Using the notation of the angled brackets 



as the thermal average, it can be shown that the zeroth moment Mo — (j ^2( r )^2(r)(iiy. 
For the first (see also Ref. [l?], Q) and the second moments we have 

Mi = - (J dr#M(rVV 2 (r)) + 5eM (3) 
M 2 = Q drg 2 Af?(r)Af 2 (r)^ + 28eM x - 5e 2 M . (4) 

Here, ftf a (r) = ^ r ^(r)^ cr (r) and 5e = (ji 2 — £2 — A*3 + £3). We restrict our study to the case 
where the chemical potentials measured with respect to the corresponding Zeeman levels 
are the same fi = fi a — e a and we have Se = 0. In this case, the atomic populations are 
balanced. The clock shift u s and the width D of the response peak are obtained from 



M 1 

oo s = — — , D 
M 



M 2 


■Mi 


M 


[m q \ 



(5) 



Thus, u> s (a s = 00) is the same as the difference of the clock shifts 5u s = j^{Mi(a s = 
00) — Mi{a s = 0)). M , u s and D can be interpreted as the parameters of the Gaussian 



it [l9[] to the response function: Sq(uj) = M /\/27rD 2 e _< - w_a;s - 1 ^ 2D \ From the experiments 
2I, this fit seems to be sufficiently accurate assuming that there are no impurities such as 
free atoms. For this case, Mi and M 2 are positive for all T. When g = (or l/a s = —00), 
the response peak has zero shift u s = and the width D = 0. uo s measures the energy 
transfer per particle in the coherent rotation of the initial 1 12} state into the |l/3) state with 
\(5) being a point in the pseudospin Bloch sphere. However, since \(3) is not an eigenstate of 
the Hamiltonian the coherence is lost due to the interactions with subsequent increasing of 
D [l?]]. Thus we have non-zero u s and D > for the interacting system (g = g 12 7^ 0) even 
when the level couplings #i 3 = g 2 3 = 0. 

Method: For the lattice model, it can be shown that around the resonance (a s ~ ±00) the 
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sign problem can be avoided with the finite and negative regulated coupling g. The complex 
phase problem returns when the molecule size (~ a s ) in the BEC regime becomes smaller 
than the lattice spacing SI and the coupling constant q becomes positive. For the finite T 
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21] . Here, the thermal operator e" 



HP 



calculations we use the grand canonical formalism 
can be applied in the single particle basis after Hubbard- Stratonovich (HS) transformation. 
The resulting multidimensional integration of the auxiliary variables is evaluated by the 
Monte Carlo methods. Typically, the inverse temperature f3 is sliced into a few hundreds of 
smaller steps 8(5. 



00 -x(r) 2 /2 

e -g J ^(r)*! (!■)*+ (r)* 2 (r)dr5/3 _ TT f ^(^f e a ! o a :(r)-a ; g/2 e -[(a : (r)-xo)-/^i+^](* t 1 (r)4ri(r)+*+(r)* 2 (r)) 

K L J a/2T 



2tt 

— oo 

(6) 

The thermalization of the stochastic samples can be optimized by shifting the center of the 
auxiliary fields. The mean shift of the auxiliary fields can be derived from the minimal 
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condition of the weight function. Thus the shift of the fields becomes x$ = ^/—g5(3n 
with n = kp/(3ir 2 ). In the practice, some freedom is given to the choice of Xq within a factor 
of order one. This shift of the auxiliary field modifies the HS transformation as shown in the 
Eq. El The Monte Carlo integration with the probability density given by Eq. [6] is carried 
out in the following steps: Firstly, the value of the field x(r) is tentatively given by direct 
Gaussian sampling. Secondly, the operator part e-K-'W-^JV^+flWaKvlWViW+vJWtfaW) 



jility density associ- 
while the factor 



is expressed in the single particle orbitals. Finally, the resulting proba 
ated with the propagator is evaluated as the determinant of a matrix 
e x x(r)-xl/2 con tributes as a simple numerical factor to the total probability density. This 
configuration is accepted or rejected by the Metropolis algorithm. In order to obtain the 
acceptance probability of 30 ~ 70%, only a small fraction of the field variables are updated 
in each trial. At T ~ O.lep just 0.5% of the field variables are updated while at T ~ 0.4ep, 
~ 3% of the auxiliary variables can be updated in each step. The benefit of shifting the 
center of the auxiliary fields is that the samples start out closer to the thermalized configu- 
ration avoiding long initial thermalization runs. However, this shift has negligible effects in 
reducing the statistical fluctuations of the thermalized samples. Nor does this shift appear 
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FIG. 1: (Color online) Here the energy per particle E in the units of Efg = 0.6cf is shown as 
circles (dotted line). All the other quantities are in the units of ep. The dashed line corresponds 



2i. 



to the chemical potential \i. Both E and \x are in good agreement with the known results 
The triangles(dot dashed line) represent the RF clock shift co s (a s = oo) and the diamonds (line) 
the width D. iv s has a 'bump' at T < 0.2eF and then remains more or less constant for higher 
temperatures. The width D of the response peak is almost featureless. It remains more or less 
constant D ~ 0.7e^ at all T. These results were obtained for 7 3 and 9 3 lattice boxes with the filling 
density of ~ 0.09. The smooth fitting curves were generated by Bezier algorithm. The size of the 
symbols correspond approximately to the statistical errors. The vertical dotted line corresponds 
to the probable position of T c 1(3) 11]. 



to ameliorate the sign problem of the unbalanced system with \i\ ^ fi2- The Monte Carlo 
simulations produce as output the one-body density matrix n CT (k, k', {x(r)}) according to 
the sampled configuration {x(r)}. Through Fourier transform, it can be connected to other 
quantities in the coordinate space. 
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251 ] . the frequency 



(7) 



Results and Discussion: In the simple BCS mean-field picture 
shift is 

,n 2 A 2 

ujmf = -.<? — ■ 

2 g'n 

where g' = -^g, n is the particle density and m* is the effective mass. In principle, g' can 
have temperature dependence through the effective mass m*(T). In Ref. [26| the calculated 
effective mass stays more or less close to ~ m in the range of temperature OTe^ ~ T < 0.9eF 
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FIG. 2: (Color online) Temperature dependence of the RF response S(uj) (by Gaussian fit described 
in the text). 
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FIG. 3: (Color online) Here the pairing gap A qp calculated using the Eq. [7] is shown for different 
effective masses: dotted line for m* jm = 1 and dashed line for m* jm = 1.43. We notice that A qp 
has a bulge at T < T c matching qualitatively with the behavior of the order parameter Ao (line 
from Ref. 9( where dynamic mean field method was used). Then for T > 0.2eF, A qp smoothly 
decays showing the transition into the gap insulator phase. A qp = 0.50(2)e^ [27] at T = is also 
shown (diamond) for the comparison purpose. 

js]. In our case, we need the ratio m*/m ~ 1.43 in order to extrapolate to the known T = 
result for the quasiparticle gap. Then we keep the same effective mass for T > 0. Our results 
are summarized in the figures [TITSl [31 and|H These results are obtained for the lattice box of 
volumes 7 3 and 9 3 with the periodic boundary conditions. Here, the finite size dependence 
is within the error bars. A qp shown in Fig. [3] is obtained from fitting u s to the Eq. [71 Even 
without considering the mass renormalization (m*/m = 1), A qp at T < T c has a qualitative 
tendency close to the known value of the A qp at zero temperature j^] . In order to optimize 
the fit, we need to adjust the mass by m* jm ~ 1.43 as earlier mentioned. We notice that 
the A qp has a bulging feature at low temperature (Pomeranchuk effect as discussed later) 
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and then decays slowly for T > 0.2ep- This is somewhat similar to the pseudogap behavior 
discussed in Ref. jf], Q]: A gp = A 2 SC + A 2 pg where A sc = A is the superconducting order 
parameter with non-zero values at T < T c and A pg is the pseudogap which is non-zero for 
T < T* with T* > T c . However, in our case there is no clear evidence of T*. In fact, at 
T > T c we enter into the regime of the gap insulator phase as discussed in Ref. The 
width of the response D shows the broadening effect due to the decoherence as discussed 
earlier (Fig. [Tf2]) . D remains more or less constant in the range of the temperature that we 
considered. In the experiment of reference [3j the population density of the K 40 atoms in 
the state 3 is measured as a function of the energy uj and the momentum k of the ejected 
particles. In this case, the momentum contribution by the RF signal is considered negligible. 
Thus the measured momentum is a good estimate of that of the non-perturbed system. Here, 
the measured A qv j [i « 0.75 at T « T c while from our estimates this ratio is ~ 1. In this 
experiment, however, the trap is spatially inhomogeneous causing the hp to vary locally. In 
another experiment (see Ref. 4J), by bimodal spectral response, the ratio A qp /fi « 1 has 
been measured at T < T c in closer agreement with our theoretical results. 

In comparison to the earlier theoretical works without the final state interaction effect 



(Ref. [6|, |28|, |29[ ) , both the RF clock shift and the width of the response function are found to 



be larger. Also, in our case u> s does not approach zero value in the studied temperature range. 



Unlike in the two channel models 



we omitted the free atom contributions 



that produce a sharp peak at the zero frequency shift. 

Connection to the Pomeranchuk effect for the lattice fermions was made in the Ref. 
[13I . [14J]. Here, anomalous double occupancy of the repulsively interacting fermions in the 
two dimensional lattice at half filling was discussed. Analogously, for the unitary limit Fermi 
gas we can relate the double occupancy to the RF clock shift by From the Fig. [TJ we 

notice an enhancement of the density overlap (double occupancy) at T < 0.2ej?. In this case, 
the adiabatic behavior of the temperature for small changes of the coupling g can be obtained 
from |jj = _2g T c (T) 9uJ q^ ■ The change in the coupling constant is 5g > (while g < 0) when 
the unitary limit is crossed from the BCS to the BEC side and an enhancement of the 
adiabatic effect that lowers the temperature occurs at T > T c (Fig. @|. This is qualitatively 
different from the 8|, [9j where the non monotonic uj s (T) can produce temperature increase. 
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FIG. 4: (Co 
in the Ref. 



or online) In the plot 1), the entropy per particle s is shown using the method described 
11]. The plot 2) shows the corresponding heat capacity calculated from c = The 



adiabatic behavior of the temperature for small changes of the coupling G = g/5l s is shown in the 
plot 3). The quantities shown here are in dimensionless units. 

Thus, precise measurement of the temperature at around T « T c with adiabatic changes 
of the resonant magnetic field can lead to the verification of the temperature dependence 
of uj s (T) (and also that of A pg (T)). Since there is no transfer of heat during this adiabatic 
process, this is not a cooling nor a hea ting effect. Similar treatment within the mean field 



formalism is also discussed in the Ref. 



3l|. 



Concluding Remarks: In summary, we have studied the RF response function of the 
unitarity Fermi gas in a fully three dimensional system by a numerical method free of any 
uncontrolled approximation. We found an unambiguous signature of the pairing gap at low 
temperatures and also the pseudogap behavior at T > T c that seems to differ qualitatively 
from some of the existing works. We described the adiabatic temperature effect which 
could lead to the experimental confirmation of the pseudogap features. We also found 
that within the quasiparticle picture a rather heavy effective mass has to be assumed. For 
the calculations, we relied on the grand canonical formalism with small statistical errors. 
There are finite temperature canonical formalisms [32] where, at least in principle, the usual 
odd-even staggering of the energy as a function of the particle numbers can be used to 
extract the pairing gap. However, the errors are known to be larger and the computational 
demands much higher. We have shown that the RF spectroscopy is a useful way to observe 
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